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ABSTRACT 

We apply the formalism of spherical harmonic decomposition to the galaxy density 
field of the IRAS PSCz redshift survey. The PSCz redshift survey has almost all- 
sky coverage and includes IRAS galaxies to a flux limit of 0.6 Jy. Using maximum 
likelihood methods to examine (to first order) the distortion of the galaxy pattern 
due to redshift coordinates, we have measured the parameter (3 = rj 6 /6. We also 
simultaneously measure (a) the undistorted amplitude of perturbations in the galaxy 
distribution when a parameterised power spectrum is assumed, or (b) the shape and 
amplitude of the real-space power spectrum if the band-power in a set of passbands 
is measured in a step-wise fashion. These methods are extensively tested on a series 
of CDM, ACDM and MDM simulations and are found to be unbiased. 

We obtain consistent results for the subset of the PSCz catalogue with flux above 
0.75 Jy, but inclusion of galaxies to the formal flux limit of the catalogue gives varia- 
tions which are larger than our internal errors. For the 0.75 Jy catalogue we find, in 
the case of a parameterised power spectrum, (3 = 0.58 ± 0.26 and the amplitude of 
the real space power measured at wavenumber k = O.l/iMpc" 1 is Ao.i = 0.42 ± 0.03. 
Freeing the shape of the power spectrum we find that (3 = 0.47 ± 0.16 (conditional 
error), and Ao.i = 0.47 ± 0.03. The shape of the real-space power spectrum is con- 
sistent with a T = 0.2 CDM-like shape parameter, but does not strongly rule out a 
number of other models. Finally by combining our estimate of the amplitude of galaxy 
clustering and the distortion parameter we find the amplitude of mass fluctuations on 
a scale k = 0.1 /i _1 Mpc is A p = 0.24il^ ' 6 , with an uncertainty of 50%. 

Key words: cosmology: large-scale structure of Universe 
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1 INTRODUCTION 

The pattern of galaxy clustering in redshift surveys is dis- 
torted because of peculiar (non-Hubble) velocities which 
contribute to the line of sight velocity measured by the 
galaxy redshift. This distortion makes the interpretation of 
the clustering more difficult, and must be corrected for when 
measuring statistical properties of the pattern, such as the 2- 
point correlation function or power spectrum. However the 



distortion itself contains important information about the 
motion of matter in the universe. Because the effect is only 
radial, the measured clustering pattern becomes anisotropic. 
The degree of anisotropy is determined by the distortion 
parameter^ [3 = fi°' 6 /6, where Q, is the cosmological mass- 

* A distinction must be made between the distortion parameters 
for different galaxy types, since their relative biases are known to 
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density parameter and b 2 = P(fc) ga i/P(fc) ma ss is the bias pa- 
rameter, relating clustering in mass to clustering in galaxies. 

By measuring the degree of distortion one can hope to 
constrain f3, and hence the density parameter and bias. In 
addition, measurement of the real-space (i.e. undistorted) 
power spectrum can put strong constraints on theories of 
structure formation. A major aim of large-scale redshift sur- 
veys is to constrain the shape and position of the break in 
the power spectrum as it turns over into the primordial spec- 
trum. This would demonstrate a strong link between the 
clustering responsible for the microwave background fluctu- 
ations and that found in redshift surveys, as well as putting 
strong constraints on the slope of the connecting spectrum. 
This in turn provides constraints on inflationary models for 
the generation of structure. 

A comprehensive overview of the literature to date con- 
cerning the analysis of redshift space distortions is given by 
Hamilton (1997). Here we shall outline some of the main 
results leading up to this analysis. 

Jackson (1972) first noticed the distortion when view- 
ing a cluster of galaxies in redshift-space, where peculiar 
velocities are generated by virial motion. Kaiser (1987) de- 
rived the main features of the redshift distortion in the linear 
regime, simplifying the analysis by assuming plane density 
waves and a distant observer approximation. The major re- 
sult was that power is boosted in redshift space by a factor 
(1 + fin\) 2 , where fik is the cosine of the angle between 
wavevector and line of sight. Hamilton (1993) suggested ex- 
panding the distorted power in Legendre polynomials, since 
the main distortion is quadrupole in nature. Taking the ratio 
of quadrupole to monopole power yields a function in terms 
of (3 only, and one could hope to measure the distortion inde- 
pendently of power. Applying this to the IRAS 2Jy survey 
Hamilton (1993) found (3 = 0.69 ± 0.26. Hamilton further 
applied this method to a merged QDOT + 1.2 Jy redshift 
survey and found /3 = 0.69 ± 0.20. However Hamilton also 
suggested that the merged survey was systematically biased 
due to a difference in clustering between the near and far 
regions of the QDOT sample. 

Cole, Fisher & Weinberg (1994) and Cole, Fisher & 
Weinberg (1995) tested the quadrupole to monopole esti- 
mator on simulations and applied it to the IRAS 1.2Jy and 
QDOT surveys and found f3 = 0.52±0.15. Peacock & Dodds 
(1994) used a radially projected Gaussian smoothing func- 
tion to model the effects of nonlinearity in the small scale 
velocity field. They found that (3 = 1.0 ±0.2 for IRAS galax- 
ies, using the APM survey to constrain the real-space power 
spectrum (Baugh & Efstathiou 1993), and allowing a rela- 
tive bias between the galaxy samples. 

Taylor & Hamilton (1996) and Fisher & Nusser (1996) 
suggested that the main contribution to nonlinearity in red- 
shift space was not from an incoherent velocity dispersion, 
but rather from coherent infall in the translinear regime. Us- 
ing the Zel'dovich approximation they calculated the scale 
dependence of the quadrupole to monopole ratio, finding 
that the amplitude was still only a function of /3. Taylor 
& Hamilton (1996) also found a dependence on the initial 
power spectrum, and the scaling behaviour of the wavenum- 



differ. Since we are concerned with IRAS selected galaxies in this 
paper, we shall reserve (3 and b to refer to IRAS galaxies. 



ber at which the quadrupole is zero. They applied these 
results to the merged QDOT + 1.2Jy surveys and found 
f3 — 0.6 ± 0.2, assuming a local spectral index of n = —1.5, 
while Fisher & Nusser (1996) found (3 = 0.6 ± 0.2 from the 
1.2Jy survey alone. 

Fisher, Scharf & Lahav (1994) and Heavens & Tay- 
lor (1995) both dropped the plane parallel approximation 
and used a spherical harmonic decomposition to match the 
spherical nature of the IRAS redshift surveys. Fisher, Scharf 
& Lahav (1994) used a Gaussian window function in the ra- 
dial direction and, applying a maximum likelihood formal- 
ism, found (3 = 0.96 ± 0.19 for the 1.2Jy survey. Heavens & 
Taylor (1995) used spherical Bessel functions to decompose 
the density field radially, since these combine the property 
of being orthonormal and eigenfunctions of the Laplacian, 
with the advantage that each mode only picks up power from 
a narrow range of wavenumbers. In addition they applied a 
radial smoothing to correct for the effects of nonlinearity. 
Applying a maximum likelihood analysis they measured si- 
multaneously (3 and the amplitude of a fixed power spectrum 
(characterised by the fractional r.m.s. in counts in 8/i _1 Mpc 
spheres), finding (3 = 1.1±0.5 and as = 0.68±0.05, where we 
quote marginal errors. Ballinger, Heavens & Taylor (1995) 
extended the analysis by dropping the fixed power spectrum 
and measuring the real-space band-power in a series of pass- 
bands in a step-wise fashion. With a free power spectrum, 
they found (3 = 1.04 with a conditional error of 0.3. 

In this paper we shall apply the methods of Heavens & 
Taylor (1995) (hereafter HT) and Ballinger, Heavens & Tay- 
lor (1995) (hereafter BHT) to measure the degree of distor- 
tion and the amplitude and shape of the undistorted power 
spectrum of the IRAS Point Source Catalogue redshift sur- 
vey (PSCz; Saunders et al. 1998). 

The contents of the paper are as follows. In Section 2 
we discuss the methods used, including some technical im- 
provements to the original analysis. In Section 3 we exten- 
sively test the methods on a suite of simulated PSCz redshift 
surveys. In Section 4 we give a brief outline of the PSCz 
catalogue itself and highlight the checks made on possible 
sources of contamination to the present analysis. In Section 
5 we present the results of the analysis for the redshift distor- 
tion and amplitude of the real-space power spectrum, and 
the shape of the real-space power spectrum measured in a 
series of passbands. In Section 6 we discuss the implications 
of these results and present our conclusions. There are also 
two technical appendices outlining the spherical harmonic 
formalism, deriving the main equations used in the analy- 
sis and detailing the construction of the mixing matrices. 
We also explain, in more detail, the technical improvements 
made to the analysis in the present work. We begin by dis- 
cussing the methods to be used. 



2 METHODS 

In this paper we present two complementary analyses of the 
redshift distortions in the PSCz survey. The first analysis, 
described in detail in Section 2.1, is a two parameter max- 
imum likelihood analysis of the data to obtain the value 
of (3 = Q, os /b and the amplitude of the real space power 
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spectrum at a wavenumber k — 0.1 /iMpc -1 ^ (assuming a 
parametrised form for the shape of the power spectrum). 
The details of this analysis follow HT with some technical 
improvements to the method concerning how the sky mask 
and Local Group motion are dealt with. These modifica- 
tions are explained in Section 2.1 and in Appendix A. The 
second analysis, described in Section 2.2 follows BHT and 
is an extension of that described above. Here we perform a 
likelihood fit for (5 and the full real space power spectrum, 
sampled in a series of passbands. This has the advantage of 
relaxing the assumption of a fixed spectral shape. We begin 
by discussing the galaxy density fluctuations in redshift- 
space. 



2.1 Density fluctuations in redshift-space 

Following HT, the density field of the galaxy distribution 
p(s) is expanded in terms of spherical harmonics, Yf m , and 
a discrete set of spherical Bessel functions, je, 



Ptmn = Ci„ Id S p(s)w(s)ji (k(„s) Yt m 



(1) 



where w(s) is an adjustable weighting function and s is the 
redshift-space position variable. The inverse transform is 



p(s) = c ln p lmn j l (k ln s)Y lm (6 , 4>), 



(2) 



where the C( n are normalization constants and ki n are dis- 
crete wavenumbers (see Appendix C of HT and Appendix 
A in this work for a definition). 

The domain of the integral in equation ([j]) can be trans- 
formed from redshift space s to real space r by continuity, 
d 3 s p(s) = d 3 r p(r), leaving only the argument of the spher- 
ical Bessel function as a function of the redshift-space dis- 
tance s. Noting that redshift distortions only affect radial 
coordinates, the transformation can be written 



s(r) — r + 



u(r) 
1h 



(3) 



where u = r.(v(r) — wlg) is the radial component of the 
peculiar velocity field in the Local Group frame and Hq is the 
Hubble parameter. Note that the expansion of the density 
field is done in this frame, since we need to assume that 
it | <C Hor, which will only be the case in the near-field 
limit if we choose a frame in which v — V{ lame — ► as r — > 0. 
In linear theory the velocity field is related to perturbations 
in the mass density by the continuity equation 



V.v = 



(4) 



where ft 6 » dlnS p /dln a is the linear growth index and S p 
is the fractional matter overdensity. Expanding the continu- 
ity equation in spherical harmonics we find that the radial 
velocity field can be expressed as 



u(r) = n 



06 E 



~j'e{k en r)Y e 



(5) 



this paper we shall use h 
is the Hubble parameter. 



Throughout 
//o/lOOkms" 1 Mpc -1 , where H 



where j[{z) = dje(z)/dz. 

Expanding equation ([j]) to first order in u — r — s we 

find 

plmn = (P0) lmn + E W "' m ( $ "" +P V U" ) <Wn', (6) 
E'm'n' 

where the mean coefficients {po)i mn are given in Appendix 
A, and include a term for the Local Group velocity, although 
the effect of this for I 7^ 1 is very small. The transition 
matrices W , 3>, and V describe the effects of the sky mask, 
the radial selection function and the first order redshift space 
distortion, respectively. The transition matrices are derived 
and defined in Appendix A. 

In HT it was assumed that the sky mask was az- 
imuthally symmetric. This simplified the analysis since the 
mask mixing matrix, W , was then real. Given the higher 
surface density of galaxies in the PSCz we drop this assump- 
tion and use the exact PSCz mask (see Figure Q). Also, we 
do not apply an artificial galactic cut, as was done conser- 
vatively for the 1.2Jy survey (HT). The mask used for this 
analysis is described briefly in Section 4 and in Appendix B 
we describe in detail the methods used for mask generation 
and application. 

Throughout this paper we have also applied a radial 
correction for the effects of small-scale nonlinear peculiar 
velocity fields (see HT, Section 4.4 and Appendix B, therein, 
and Appendix A in this work) . This assumes that the effects 
of the 'fingers of God' due to virialised clusters and redshift 
measurement errors can be modelled by adding an incoher- 
ent random displacement to each galaxy radial distance 



s' = s + e(r), 



(7) 



where (e) = 0, (e(r)s(r')) = {a 2 v /iH^)8 D {r - r') and a v 
is the three-dimensional velocity dispersion. If the distribu- 
tion of random displacement is assumed to be drawn from 
a Maxwellian distribution we should ensemble average over 
this distribution when forming the correlators of the density 
field, {S5'). But since the random displacements are incoher- 
ent we need only calculate the effects on the mean value of 
modes, since (XY) = (X)(Y) if X and Y are uncorrelated. 
The averaging effect on a given mode can be expressed as 
a convolution of the density field by a scattering matrix, S 
(see Appendix A for a definition) . 

The statistical properties of the redshift space density 
field can be encapsulated by the mode-mode correlation 
function 

(DD ] ) = i ^ SW($ + (3V)P($ + pVfW^S 1 + N (8) 

where the observable data vector, D, is given by D( mn = 
[pimn — (po)tmn]/ P and p is the mean number density, and 
N is the shot-noise contribution to the covariance matrix. 
We have assumed that the underlying Fourier modes of the 
density field are statistically isotropic and homogeneous: 



(<y(fc)<T(fc')> = {2^fP{k)8 D {k - k'), 



(9) 



where 8d is the Dirac delta function. Since the shot-noise 
term contains a part linear in the matrix W . the mode- 
mode correlator is complex for a non-azimuthally-symmetric 
mask. This can be easily overcome by decomposing the 
transformed density field into real and imaginary parts: 



D = $hD + iQmD. 



(10) 
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This is convenient since our hypothesis that D is a Gaus- 
sian random variable requires that the terms $hD and ^SmD 
are also Gaussian random variables. In the limit of full sky 
coverage these terms would also be uncorrelated, but the 
presence of a non-azimuthally symmetric mask mixes real 
and imaginary parts. Even so the likelihood function for the 
observed redshift density field can be written in Gaussian 
form: 



£[D\p,P(k)] = 



1 



r exp 



2 



C X D 



(11) 



(2n) N / 2 \C\ 1 / 2 

where D = ($hD ,$smD) is a data vector of real and imagi- 
nary parts of the density modes, with dimension its being set 
by the upper limit to the wavenumber. The data covariance 
matrix is 



(12) 



where CW = (HkDSkD*), Csfessn = CUhDSmD*} and so 
on. From the symmetry relations given in Appendix B we 
see that this is a real, symmetric matrix. 

Another refinement to the analysis of HT is the inclu- 
sion of m = modes. Since these are real modes, which 
have twice the power as the m 7^ modes, we treat these 
separately. 

Finally, we have also included a correction for the effect 
of the Local Group motion. In HT and BHT, we simply 
removed the dipole contribution to the harmonics. While 
this is exact for full sky coverage, for incomplete sky coverage 
there is a mixing of modes and the dipole leaks into higher 
modes, although at a level which was negligible for the PSCz 
mask used. Here we account for that leakage. Again technical 
details can be found in Appendix A. 

The free parameters for the likelihood are then the red- 
shift distortion parameter (3 and the real-space power spec- 
trum P(k). In this paper we shall consider simultaneous 
measurement of 

(a) j3 and the normalisation of P(k) for a fixed spectral 
shape, 

(b) j3 and the shape and amplitude of P(k). 

The first method maximises the likelihood over a 2- 
dimensional parameter space, while the latter requires a 
maximisation over a 7-dimensional parameter space - one 
dimension for j3 and 6 dimensions for the step-wise power 
spectrum sampled uniformly in logarithmic pass-bands. In 
Section 2.2 we shall discuss the step-wise maximum likeli- 
hood method for measuring the real power spectrum. Before 
doing so there are two issues we shall deal with here - that of 
the correct weighting function to use in the transformation 
equation (|l|) and the choice of how to measure the normali- 
sation of the power spectrum. 

The minimum variance weighting scheme to use in the 
transformation equation (|l]) has been discussed by Feldman, 
Kaiser & Peacock (1994), HT, and Hamilton (1997). All 
three have found that for an undistorted, flux limited red- 
shift survey the optimal weighting function for measuring 
power spectra is 

W{r) = 1 + (13) 
where <j>(r) is the survey selection function in real space and 
P w {k) is the expected power spectrum. As can be seen one 



has to make an initial guess for the power spectrum to insert 
into the weighting function. This is not important: a poor 
choice would lead to larger than necessary error bars, but 
should not bias the solution. In any case, if one was con- 
cerned, one could iterate. We have chosen to use af = 0.2 
power spectrum for the analysis of the PSCz, where Y is the 
CDM-like spectral shape parameter. In the tests on mock 
PSCz catalogues (Section 3) we used the P w (k) appropri- 
ate to the simulations. In both cases we have used the full 
fc-dependence in the weighting scheme (cf HT, who used a k- 
independent P(k) in the analysis of the 1.2Jy survey). Note 
that this weighting function is strictly only optimal for es- 
timating power spectra, and when (3 = 0. HT derived an 
approximate weighting scheme optimised for measuring /3, 
with the main difference being to correct for the misplace- 
ment of weights in redshift space and to boost the power 
by a factor (1 + 0). Since these seem to have little effect 
on the maximum likelihood solution, we used the weighting 
function optimised for power. 

In the two-parameter fits, we have chosen to express the 
power normalisation in terms of Ao.i, where 



Al = P(k)k 3 /(2n 2 ) 



(14) 



is the contribution to the overdensity variance per ln(fc). 
This choice is motivated by the fact that we choose an upper 
wavenumber limit k ~ O.l/iMpc -1 to avoid nonlinear scales, 
and most of the modes are concentrated near the upper 
wavenumber limit. Hence it is the power at this wavenumber 
which is most strongly constrained and is the most insensi- 
tive to the shape of the power spectrum assumed. Since it 
is common to normalise spectra by the variance in spheres 
of radius 8 /i _1 Mpc, 



2 

cr 8 



d ±A\k) 



3ji(fcr 8 ) 
kr$ 



(15) 



we have fitted the relationship between erg and Ao.i as a 
function of the CDM shape parameter, F. A fit accurate to 
a few percent is given by 



ffg = 0.616 



0.84 1 



r 

02 



Ao.i 
0.5 



(16) 



A further consideration is the validity of the linear bias 
model used to relate the galaxy distribution to the peculiar 
velocity field, which appears in /3. While a wavenumber of 
k = 0.1/iMpc _1 corresponds to a scale of about 60 /i _1 Mpc, 
and one would hope biasing could be described by a local 
and linear model, we have no concrete guarantee that bias- 
ing is not still scale-dependent or non-local, as in the case 
of cooperative biasing (Bower et al. 1993). However realistic 
models of biasing mechanisms do seem to give some sup- 
port to local biasing that is approximately linear on these 
scales (Kauffman, Nusser & Steinmetz 1997, Mann, Peacock 
& Heavens 1998, Heavens, Matarrese & Verde 1998). 

Having described the main formalism for the measure- 
ment of distortions and real space power in this section, we 
now turn to the methods used to measure the shape and 
amplitude of the power spectrum in a stepwise fashion. 

2.2 Stepwise maximum likelihood 

This method fits the shape of a function in an essentially 
model-independent way - see Efstathiou, Ellis & Peterson 
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(1988). The function is divided into bins, and the amplitude 
of these bins is allowed to vary. 

We exploit the fact that we have a good estimate for the 
power spectrum from previous studies, Po(k), so we write 

P(k) = P (fc) x 10*, (17) 

where the original power spectrum (CDM with V — 0.2 and 
Aq.i = 0.18) is multiplied by the free function 10^. The 
subscript i = 1, 2, 3, . . . labels each wavenumber bin corre- 
sponding to the range fcf m < k < fcf lax . The parameters 
are introduced in this exponential form to ensure that the 
measured power spectrum is positive. 

The parameters for the maximum likelihood fit are the 
distortion parameter (3 and the coefficients i = 1, . . . , n. 
This construction essentially amounts to specifying the 
power spectrum at a series of wavenumbers, but the form 
chosen allows us to start the maximum likelihood search 
conveniently at fi = 0, since we have a fairly good idea of 
what to expect. Note that, although the form suggests that 
the break scale is fixed, adjustment of the /; can alter this 
in an essentially arbitrary way, which is the main reason 
for parametrising the power spectrum in bins. In HT, and 
in Fisher, Scharf & Lahav (1994), the power spectrum was 
parametrised to have a specific shape, whereas here it has a 
free form (apart from bin discreteness). 

We begin by expanding the covariance matrix of equa- 
tion @ in powers of /3; 



C = J2 l 3 "^ + N 



(18) 



The covariance matrix can now be constructed with lit- 
tle computational cost as each sub-matrix can be written 
as a linear combination of contributions from each power- 
spectrum bin; 

A = 10 /l Aoi + 10 /2 A 02 + ..., 
Ai = 10 /l A 11 + 10 /2 Ai2 + ..., 

A 2 = 10 A A 2 i + W h A22 + ■ ■ ■ , (19) 

where A01 is the contribution to Ao from the first bin, and 
so on. The matrices Ao%, An and A21 need only be con- 
structed once and then the full A-matrices can be created 
each time from a simple linear sum. This requires far less 
computer time than fitting for a standard power spectrum 
parameter such as T. That would require the A-matrices 
to be constructed afresh from window function convolutions 
each time a parameter changed. 

The maximizing procedure is straightforward, the (n + 
1)D parameter space is searched until a maximum of the 
likelihood function C\D\f3, fi\i = l,...,n] is found. This 
generalizes the method of HT, where, in effect, the fi were 
all constrained to move up or down together. 

If n is high, exploring the likelihood space with a grid 
can become computationally expensive. However, the maxi- 
mum likelihood solution is easily found, since we have a good 
starting guess for the parameters, and can use a standard 
library routine to maximise the likelihood function with re- 
spect to the parameters. Ideally, we would like to construct 
complete hypersurfaces of constant likelihood, but this is 
impractical because of the dimensionality of the problem. 
Errors can be obtained by searching along the parameter 
axes until the likelihood drops by y/e, but note that this ig- 



nores the correlations between parameters. This means that 
the errors are conditional: for each parameter they are calcu- 
lated assuming that the other parameter values are known 
exactly - see Press et al. (1992; §15) for a discussion. Note 
in particular that there is a strong covariance between the 
amplitude of the power spectrum and 0, and these errors 
do not properly reflect this (see HT, Fisher, Scharf & Lahav 
1994). 

By parametrising in bins, the fitted power spectrum is 
discontinuous. In principle there might be a case for adjust- 
ing the shape of the power spectrum so that it passes contin- 
uously through the maximum likelihood points; in the event 
there is little to be achieved by this, since the error bars are 
large enough so that a smooth curve is a perfectly acceptable 
fit. 



3 TESTS ON SIMULATED SURVEYS 

The methods described in Sections 2.1 and 2.2 have been 
applied to mock PSCz redshift surveys drawn from large 
N-body simulations. Before presenting the results from the 
PSCz redshift survey, we describe the result of applying 
these methods to the mock PSCz data sets. 



3.1 Numerical simulations 

The dissipationless N-body simulations that we have used 
to create the mock data-sets consist of 3 ensembles of 9 
realizations each. The first two ensembles are of CDM-like 
universes and contain 160 3 particles within a periodic com- 
putational box of length it = 600 /i _1 Mpc. These simu- 
lations are similar to the ones used and discussed in de- 
tail by Croft & Efstathiou (1994a) but employ more par- 
ticles within a larger computation box (the simulations of 
Croft and Efstathiou use 100 3 particles within a box of 
size £b — 300 7i _1 Mpc). The simulations were run with the 
particle-particle-particle-mesh (P 3 M) code described by Ef- 
stathiou et al. (1985) and model the gravitational clustering 
of dark matter in a CDM dominated universe with scale in- 
variant initial density fluctuations. These two ensembles are 
as follows: the standard CDM model (Davis et al. 1985), i.e. 
a spatially flat universe with Qo — 1 and h = 0.5 (the SCDM 
ensemble); a spatially flat low density CDM universe with 
fio = 0.2, h — 1.0, and a cosmological constant contribution 
tt A = A/(3#o) = (1-fio) = 0.8 (the ACDM ensemble). 
The initial power spectra of the models are generated from 
the fitting function 



P(k) oc 



k 



1+ ak+ (6fe)2 + (ckf 



(20) 



where v = 1.13, a = (6.4/r) h^Mpc, b = (3.0/r) /i _1 Mpc 
and c = (1.7/r) /i~ 1 Mpc. Equation (18) is a good approxi- 
mation to the linear power spectrum of scale-invariant CDM 
models with low baryon density, fij, <C ilo (Bond & Efs- 
tathiou 1984). The parameter T in equation (18) is equal 
to 0.5 for the SCDM ensemble and F = 0.2 for the ACDM 
ensemble. For pure CDM-type models the shape parameter 
is given by F = tth. But if the model contains baryonic mat- 
ter the shape parameter can be rescaled so F = D.he~ nB ^ 2 , 
where Qb is the baryon density parameter. (Peacock 1992). 
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The third ensemble of 9 simulations is a spatially flat 
mixed dark matter model (the MDM ensemble) in which 
CDM contributes Ocdm = 0.6, baryons contribute Qb = 
0.1 and massive neutrinos contribute Q. v — 0.3. This set of 
simulations was run in a box of length 300 /i _1 Mpc with 
100 3 particles. The initial conditions for these simulations 
are generated from the power spectrum given in equation 
(1) of Klypin et al. (1993). In the MDM models, the thermal 
motions of the neutrinos are ignored and so they follow the 
evolution of a collisionless cold component with Qo = 1. 
The MDM simulations are exactly those used by Croft & 
Efstathiou (1994b). 

In the analysis presented below, we use data and 
simulated data out to a comoving distance of at least 
240 /i _1 Mpc, this means that the MDM models are not large 
enough to properly simulate the survey without a signifi- 
cant amount of wrapping of the simulation box. The larger 
CDM-like models however are superior in this respect and 
no wrapping is required to simulate the data. 

The final output times of the models are chosen to ap- 
proximately match the microwave background anisotropics 
measured in the first year COBE maps (Smoot et al. 1992) 
ignoring any contribution from gravitational waves; 



A(fc) = 180fc 2 T(k). 



(21) 



Thus the rms mass fluctuations in spheres of radius 
8 h~ 1 Mpc are as — 1 for the CDM-like models and as = 
0.67 for the MDM model. The amplitude of the tempera- 
ture anisotropics measured from the four year COBE maps 
(Wright et al. 1996) is slightly higher than from the first year 
results. However, for the purposes of testing the techniques 
described above, this discrepancy is unimportant. 



3.2 Construction of mock PSCz surveys 

We draw mock PSCz redshift surveys from the N-body sim- 
ulation distribution as follows: 

(i) An observer, situated at an arbitrary position in the 
simulation box, carves out a sphere of radius r max = 
240 /i _1 Mpc. 

(ii) Mass points are selected according the PSCz selection 
function: 



0(A) = fa 



10 



(l-a)A 



(1 + 10t a )~ 



(22) 



where A = log 10 (r/ro) and the values of the parame- 
ters are fa = 0.006794 /i 3 Mpc" 3 , r = 1.95113 /i _1 Mpc, 
a = 1.88654, rj = 4.38595 and 7 = 1.54554. These param- 
eters were derived using the methods described by Mann, 
Saunders & Taylor (1996). 

(iii) Mass points are moved to their redshift space posi- 
tions, s, via the equation 



1 + 



(v - vq) 
H r 



(23) 



where we have full information on the velocity v of each 
particle. The observer velocity vo is taken to be zero for the 
simulations, but is more correctly the Local Group veloc- 
ity in the survey (so the second term in brackets does not 
formally diverge as r — » 0). 



(iv) Cartesian x,y, z coordinates are converted to galactic 
latitude and longitude and a PSCz sky mask (see Section 
4.1) is applied to the distribution. 

(v) Mass points are drawn at random to have the correct 
selection function (p^). 

(vi) Mass points in the masked region of sky are dis- 
carded. The mask used here is in fact slightly less conserva- 
tive than that finally used for the real data. 

Using this procedure we obtain mock PSCz surveys with 
approximately the same density (and therefore the same 
number) of mass points as there are galaxies in the PSCz 
survey to a proper distance of r max /i~ 1 Mpc. Mass points 
are taken to represent galaxies throughout the analysis, so 
the linear bias parameter is unity. We have made no attempt 
to choose 'galaxies' from the simulation density distribution 
as this introduces extra complications and is unnecessary 
for the statistical testing purposes for which we employ the 
mock catalogues. 



3.3 Results of analysing mock PSCz catalogues 

Figures [l], [2] and |^ show the result of applying the analysis 
described in Section 2.1 to the mock PSCz catalogues. In all 
cases, we have used r max = 240 /i _1 Mpc and have analysed 
modes up to k = 0.1 /iMpc" 1 . The maximum t and n values 
are 21 and 8. For the convolutions we have included modes 
up to £ = 30, n = 20. 

For the calculation of the scattering matrix (equation 
I52} Appendix A), the three-dimensional velocity dispersion 
a v was measured from the simulations for each model. We 
have used the values: a v — 1160 kms^ 1 , 587 kms" 1 and 
859 kms" 1 for the SCDM, ACDM and MDM models respec- 
tively. 

In the analysis of the CDM-like mock catalogues we 
have assumed a real-space power spectrum of the form given 
in equation (^o|) . For the MDM model, a power spectrum of 
the form of equation (1) of Klypin et al. (1993) was used. 
These power spectra are needed for the calculation of the 
covariance matrix (equations 52 - 54, Appendix A) as well 
as for the weighting function (equation Il3f) . 

Figures [jj ^ and ^| show contours of likelihood in the 
- A0.1 plane for the SCDM, ACDM and MDM mock cat- 
alogues respectively. The contours are plotted at intervals of 
0.5 in log-likelihood and the :r-axis is labelled at intervals 
of 0.1 in A0.1. The true values of /3 and A0.1 are indicated 
by points. Figure ^ shows the ensemble average of the indi- 
vidual likelihoods and demonstrates the correct recovery of 
the parameters 

The average recovered value of /3 and A0.1 over the nine 
realizations of the PSCz survey for each of the cosmological 
models, together with the true values of these parameters 
are presented in Table [j]. The error on the recovered value 
of j3 and A0.1 is the error on the mean value (obtained from 
averaging the estimates, rather than from the ensemble av- 
eraged likelihood plots of Fig ^) multiplied by y/9 i.e. this 
is the error appropriate to the recovered values of (3 and 
A0.1 from a single realization of the PSCz survey. This er- 
ror bar provides us with an 'external' estimate of the error. 
The formal 'internal' errors on the values of [3 and A0.1 are 
given by the projection of SlnC = —0.5 onto the parameter 
axes, i.e. the projection of the first plotted contour. Evalu- 
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Figure 1. Figure shows contours of likelihood in the j3 — Ao.i plane for 9 realizations of the PSCz redshift survey drawn from SCDM 
N-body simulations. Points show the true parameter values. The contours are plotted at intervals of <51n£ = 0.5. The joint distribution 
has (/3, Ao.i) within 2.3 contours with 68% probability, and within 6.17 contours with 95% probability. The 1-er error bars on the 
parameters individually are obtained by projecting the inner contour onto the axes (Press et al. 1992). 



ating these points for each of the plots shown in Figures |lj 
| and | gives the internal error estimates for j3 and Ao.i, 
together with the one sigma range on this error. These er- 
rors are shown in columns 3 and 6 of Table |l| for (3 and 
Ao.i, respectively. Note that the parameter estimates for 
the external error determination were made from likelihood 
surfaces which were rather crudely sampled, so the mean of 
these estimates is artificially close to the true solution. 

The recovered values of (3 and Ao.i are in very good 
agreement with the true values. In each model, the correct 
values of (5 and Ao.i are recovered to within the la error. 
The internal error measured from the projection of the first 
contour level provides a reasonable estimate of the true er- 
ror, as can be seen from the fact that the external and in- 
ternal error estimates are in good agreement. The only case 
where this is not true is for Ao.i in the MDM model where 
the external error appears to be too large. This is possibly 
due to the fact that these simulations are not large enough 



to model the PSCz redshift survey out to a proper distance 
of r max = 240 ft _1 Mpc. 

It is remarkable that the correction applied here for the 
non-linear effects introduced by the small-scale velocity field 
is adequate even for the SCDM model, which has a very high 
three-dimensional velocity dispersion (a = 1160 kms -1 ). 
Without this correction we would expect the recovered value 
of f3 to be significantly lower than the true value as non- 
linear ("finger of god") effects become important even for 
relatively small wavenumbers (k ~ 0.06 /iMpc -1 - see e.g. 
Tadros & Efstathiou 1996). 

We have also tested the stepwise maximum likelihood 
technique using the mock PSCz surveys described above. 
Figure ^| shows the results for individual SCDM mock cata- 
logues and Figure ^ shows the ensemble averaged results for 
all of the SCDM, ACDM and MDM models. As expected in 
individual simulations the uncertainty grows with decreas- 
ing wavenumber, but the ensemble averaged results show 
that the stepwise maximum likelihood method is unbiased. 
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Figure 2. Contours of likelihood in the /3 - Ao.i plane for 9 realizations of the PSCz redshift survey drawn from ACDM N-body 
simulations. Again, as in all the following figures, contours arc plotted at intervals of <51n£ = 0.5. Points show the true parameter values. 



Note in particular that the longest measured wavelength 
bin shows no significant bias, since we do not include the 
t = mode. This has caused problems in the past for power 
spectrum analysis where the mean of the survey is estimated 
from the survey itself, suppressing the longest wavelength 
modes (Tadros & Efstathiou 1996). 

From these tests on simulated PSCz redshift surveys 
we can be confident that the analysis described in this pa- 
per provides an unbiased estimate of the parameters j3 and 
Ao.i. Furthermore, the errors obtained by examining the 
likelihood contours provide a reasonable estimate of the true 
error on j3 and Ao.i. 



4 THE IRAS PSCZ REDSHIFT SURVEY 

We are now in a position to apply the analysis described 
above to the newly completed IRAS PSCz redshift survey. 
In Section 4.1 we briefly describe the survey itself, and in 
Sections 4.2 & 4.3 we test for various systematic effects in 



the catalogue that could bias the recovered value of f3. The 
results of our analysis are presented in Section 5. 

4.1 The Dataset 

The IRAS PSCz survey is the largest and deepest all sky red- 
shift survey of galaxies to date. It was conducted as a follow 
up to the 1 in 6 sparsely sampled QDOT redshift survey 
(Lawrence et al. 1998). The aim was to determine redshifts 
for all galaxies in the IRAS Point Source Catalogue (PSC 
- Joint IRAS Science Working Group 1988) to a 60/im flux 
limit of 0.6Jy in the area of sky where reliable optical iden- 
tifications and redshift acquisition were feasible. An outline 
of survey construction is given in Saunders et al. (1995). 
About 4600 new redshifts were measured in the course of 
the project, the rest being taken from the literature or ob- 
tained via private communications. Details of the data re- 
duction process will be published elsewhere - (Keeble et al. 
in preparation). Figure ^| shows the N(z) distribution for 
galaxies in the PSCz redshift survey. The median redshift of 
galaxies in the survey is ~ 0.03, but there is a tail of galax- 



© 0000 RAS, MNRAS 000, 000-000 



Spherical Harmonic Analysis of the PSCz Galaxy Catalogue 9 





3 


: 




2 


5 r 




2 


: 


d 

a 


' 


5 7 


II 


1 


r 







5 : 



0.30 0.40 





3 


: 




2 


5 r 




2 


: 


d 

a 


' 


5 7 


II 




o - 







5 7 





3 


: 




2 


5 7 




2 


r 


d 
C 


1 


5 : 


II 


1 


o - 







5 7 




0.30 0.40 





0.30 0.40 0.50 0.60 




0.30 0.40 0.50 0.60 



0.30 0.40 0.50 0.60 

A„, 



Figure 3. Figure shows contours of likelihood in the /3 - Ao.i plane for 9 realizations of the PSCz redshift survey drawn from MDM 
N-body simulations. 



Table I 



Model 




2. Err or i 




4.A rec 


5. Error j 


6. Atrue 


SCDM 


1.01 ± 0.24 


0.17 ± 0.05 


1.00 


0.53 ± 0.03 


0.04 ± 0.01 


0.54 


ACDM 


0.36 ± 0.12 


0.15 ± 0.04 


0.38 


0.66 ± 0.03 


0.05 ± 0.01 


0.64 


MDM 


1.00 ± 0.60 


0.20 ± 0.05 


1.00 


0.51 ± 0.06 


0.03 ± 0.01 


0.49 



Table 1. The results of applying the analysis described in Section 2.1 to mock PSCz redshift surveys drawn from three cosmological 
models. Columns 1 and 4 give the recovered values of j3 and Ao.i, averaged over nine realizations, with errors derived from the scatter 
in the recovered values (external errors). These are to be compared with the true values of the parameters, given in columns 3 and 6. 
Columns 2 & 5 show the mean and scatter of the 1-cr projected error bars (internal errors). 



ies out to redshifts as large as ~ 0.2. The solid curve shows 
the predicted N(z) using the selection function of equation 
(p2|). The total error on the measured velocities is ~ 100 - 
150 kms^ 1 , which is negligible for this analysis. 

The PSCz was undertaken for various analyses, some 
of which demand the largest sky coverage possible. The de- 
fault PSCz mask defining the survey consists of the IRAS 
coverage gaps, the LMC and SMC, and areas estimated on 
the basis of their WOfim sky brightness (Rowan- Robinson 



et al. 1986) to have A B > 2 m , altogether 16% of the sky 
(Saunders et al. 1998). At low latitudes there is inevitable 
degradation in the survey due to incomplete identifications 
and the greater difficulty of getting redshifts. For an analysis 
such as this, where we are looking for weak clustering over 
large scales, uniformity is paramount. We therefore created 
more conservative masks, corresponding to optical extinc- 
tions of A B = 0.5 m ,0.75 m , l m , for use in this analysis. The 
results presented are all for Ab = 0.75 m (excluding 35% of 



© 0000 RAS, MNRAS 000, 000-000 



10 H. Tadros et al. 



- 




1.22±0.22 

... "\ \ 






1.18±0.23 






f-0.81±0.17 

+-'"" 






A" 
















t. 






T 








T 






31 O. 


32 O.OS O.I 0.2 ~© 


31 O. 


02 O.OS O.I 
k/h Mpc -1 


0.2 ^O. 


31 O. 


□2 O.OS O.I O. 
!</!-> Mpc"' 



* 9 



j9— 1 .37±0.23 



t. 



J 



l.OI 0.02 O.OS O.I O.S 

fc/h Mpc -1 



%01 0.02 O.OS O.I O. 

It/h Mpc - ' 



■ST s 



(?-0.87±0.20 



>.01 0.02 O.OS O.I 

k/h Mpc -1 





f-0.70±0.18 | 




j»— 1.06! 


t0.20 
















A 








f 1 
















T 










t 


l ' 





■O.01 0.02 O.OS O-l 



■O.OI 0.02 O-OS (J. 1 O.Z 

KXn Mpo-' 



■O.OI O.OZ O.OS O.I 0.2 

K/h Mpc 1 



Figure 5. The real-space power spectrum obtained from an ensemble of SCDM mock PSCz simulations (data points), with the true 
power spectrum dotted. For this analysis, r ma x = 300/i -1 Mpc and fc m ax = 0.87fe Mpc -1 . 



the sky), but are robust to which of these three latter masks 
is used. 

There is also incompleteness at high redshift, in that 
galaxies with than A^j > 19 m were not systematically pur- 
sued spectroscopically. Based on the joint optical/60/im lu- 
minosity functions of Saunders et al. (1990) and the ULIRG 
survey of Clements et al. (1996), we believe that a negligible 
number of these faint galaxies are closer than z = 0.1, for 
areas with Ab < l m . 

Within the volume defined by these cuts in area and 
distance, there are potentially variations in the IRAS data 
quality due to: the number of scans made by the satellite 
(HCONs), the South Atlantic Anomaly, Malmquist effects, 
hysteresis, confusion and noise lagging (Joint IRAS Science 
Working Group 1988). The most serious of these is the first, 
causing known incompleteness at faint fluxes. We restored 
most of this by reclaiming galaxies from the PSC reject file 
(Saunders et al. 1998). We investigated the residual qual- 
ity variations both theoretically and empirically (through 
source counts), and were unable to find evidence for any 
non-uniformities above a level of 2% in flux or 3% in den- 
sity. Several effects are expected to come in below this level, 
and so these numbers represent our overall best-guess esti- 
mates of the PSCz uniformity. Most of the concerns about 
the uniformity of the PSC affect sources near our flux limit 



of 0.6Jy, and so our final quoted results refer to a sample 
cut conservatively at 0.75Jy. 

The distribution of the PSCz galaxies used in this anal- 
ysis on the sky is shown in Figure R. 



4.2 Internal consistency checks 

The ability to make high-precision estimates of the power 
spectrum with this method allows us to look for systematic 
effects in the catalogue. We do this by making a number of 
changes in the sample and the modes analysed and looking 
for inconsistencies in the recovered power spectrum and (3. 
In this section, we give details of the various tests for sys- 
tematics we have applied to the PSCz, and in Section 4.3 
we describe the results. 

We have tested the catalogue (and the method) for sen- 
sitivity to the following: 



• Changes in the upper wavenumber limit, fc max 

• Changes to the outer boundary, r max 

• Evidence for excess angular variance (explained below) 

• Changes to the flux limit 

• Changes from point source fluxes to addscan fluxes 

• Changes between a high flux sample and low flux sam- 



ple. 
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Figure 7. The sky-distribution of galaxies in the PSCz survey with GOfim flux above 0.75Jy. Here we have applied the conservative 
sky-mask used in our analysis. No cut in redshift has been imposed. 
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Figure 4. Parameter likelihoods for ensemble- averaged mock 
PSCz redshift surveys drawn from SCDM simulations (top right), 
ACDM simulations (bottom right) and MDM simulations (top 
left). True parameter values are marked with a point. 



Figure 6. The ensemble average of real-space power spectrum es- 
timates from the masked SCDM, ACDM and MDM mock surveys, 
along with the true power spectrum in each case. The recovered 
value of (3, whose true values are 1.0, 0.38, and 1.0 are also shown. 



By 'addscan' fluxes, we mean fluxes determined from 
the coadded raw IRAS data (kindly provided by IPAC), us- 
ing software provide by Amos Yahil. The last test mentioned 
above is an important one, and needs further explanation. 
We first divide the PSCz catalogue into thin radial shells 
(width 100 km s _1 ). The median flux of galaxies in a shell 
is calculated and galaxies with flux above this value are put 
into a 'bright' catalogue (the high-split catalogue) and the 
galaxies with flux below the median go into the low-split 
catalogue. This effectively divides the PSCz by luminosity 
while ensuring that the two sub-samples have the same ra- 
dial distribution. Consequently, any inconsistencies between 
the results from the high and low-split catalogues cannot 
be due simply to changes in selection function or, equiv- 



alently, the radial distribution of points. We should point 
out that demonstrating consistency is not always straightfor- 
ward: even disjoint samples of galaxies are not independent, 
as they are (by assumption) drawn from the same underly- 
ing density field. They do, however, have independent shot 
noise. 

4.3 Results of consistency tests 

For speed, we generally chose a relatively small wavenumber 
cutoff, fc max = 0.08/i Mpc~ , when performing these consis- 
tency tests. We find no evidence for any systematic changes 
in the results for changes in r max , or fe max , or from chang- 
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Figure 8. The redshift distribution for the galaxies in the 
>0.75Jy PSCz sample, in the unmasked regions of sky. For the 
analysis, the data are also cut at various values of z max around 
0.08. The line shows the selection function used. 

ing point source fluxes to addscan fluxes. However when we 
changed the flux limit of the PSCz sample from 0.6Jy to 
0.75Jy the results for Ao.i from these catalogues (hereafter 
PSCz0.6 and PSCz0.75) differed significantly. The high/low- 
split catalogue results show this inconsistency more clearly 
as explained in the next section. 

We can make high and low-split cat alog ues from both 
PSCz0.6 and PSCz0.75. Figures | and [u] show parame- 
ter estimates for the high-split and low-split catalogues for 
these two flux limits. Figure |^ shows clear problems for 
the PSCz0.6 sample, as evidenced by the inconsistent pa- 
rameter estimation from the high and low split subsam- 
ples. PSCz0.75 on the other hand, gives consistent results. 
We have performed a counts-in-cells analysis (c.f Efstathiou 
1995), and similarly find an excess variance in the low- flux 
data. 

The PSC has known imperfections at our flux limit of 
0.6Jy. These include incompleteness in 2HCON large 
and crude 'faint flux overestimation' correction factor, and 
low and variable signal-to-noise (all Beichman et al. 1988). 
These problems have been investigated and, where possible 
corrected for in the galaxy catalogue used for this survey 
(Saunders et al. 1998). Nevertheless, the completeness and 
uniformity, which are paramount to this analysis, are certain 
to deteriorate towards the flux limit. 

Because we find that our results are sensitive to the 
flux limit we use we have conservatively chosen to limit the 
catalogue to 0.75Jy in our quoted results. 

However, despite extensive tests, we have not been able 
to find any 'smoking gun' in the sense of identifying where 
the excess clustering is coming from at lower fluxes. In fact, 
PSCz0.6 gives a result for f3 which is consistent with the 
PSCz0.75 sample; it is the amplitude which is the problem, 
being enhanced over PSCz0.75 by 15%. Given the smallness 
of the errors in the amplitude such a discrepancy is signifi- 
cant. 

We have also searched for, and failed to find, a sig- 
nature of direction-dependent flux errors, as proposed by 



Hamilton (1995). This arose from his analysis of the redshift- 
space correlation function of various IRAS surveys. He found 
an elongation along the line-of-sight on larger scales than 
plausibly accounted for by fingers-of-god. This was partic- 
ularly marked for the QDOT survey in the region beyond 
80 /i _1 Mpc, and one explanation proposed was a direction- 
dependent flux error of ~ O.lJy, varying on angular scales 
of 7°. 

We can check partly for this effect by inspecting the 
n = 1 modes for signs of excess angular variance. These 
modes, which oscillate once radially, should be enhanced, 
as flux errors would affect the distant part of the catalogue 
more than the nearby region. We see no such signature: the 
power in these modes is consistent with the model solution 
and shot noise. It should be stressed however that our large- 
scale analysis does not probe angular scales as small as 7° 
(£ = 50); it is only sensitive to variations on rather larger 
angular scales. 

We also rule out a straightforward radially-dependent 
effect as the cause of the discrepancy between the high 
and low luminosity sub-samples of PSCz0.6, since these 
sub-samples have the same mean radial distributions by 
construction. Sutherland et al. (1998) have analysed the 
redshift-space power spectrum and find the same behaviour 
as seen here, i.e. a higher amplitude when all galaxies down 
to the 0.6 Jy flux limit are included. The effect is smaller in 
redshift space, but this would be expected if the lower flux 
galaxies had an excess variance, for whatever reason. 

If it is a real luminosity effect (i.e. faint galaxies cluster 
more strongly than bright galaxies, introducing a differen- 
tial bias), it must be very large, since in the analysis, the 
samples with different flux limits differ little in the weight 
given to galaxies of different luminosity. Moreover, if it is 
a true luminosity effect we would also expect to see it in 
PSCz0.75. 

In our desire to leave no stone unturned, we also in- 
vestigated the selection function itself, as the method of de- 
termining it (Mann, Saunders & Taylor 1996) provides the 
selection function in redshift space, not real space. The ef- 
fects of this can be accounted for through a Taylor expansion 
and the addition of an extra term in the equations, but in- 
clusion of this extra term makes an infinitesimal difference 
to the parameter estimation, so we neglect it. 



5 THE REDSHIFT DISTORTION OF THE 
PSCZ AND THE AMPLITUDE OF THE 
REAL SPACE POWER SPECTRUM 

We have analysed the PSCz0.75 survey using the method 
described in Section 2.1. Here we present results for a bound- 
ary at r max = 200 /i _1 Mpc. With the sky mask as defined in 
Section 4.1 we have 7042 galaxies in the data-set to the flux 
limit of 0.75Jy. We use a wavenumber limit of fc max = 0.13h 
Mpc -1 . These modes have I up to 25, and n up to 9, but 
not all of the modes with these limits are used, as many 
fall outside the wavenumber range. We include modes up to 
I = 30 and n = 20 for the convolutions. 

For the calculation of the covariance matrices we have 
assumed a CDM-like real space power spectrum with F — 
0.2 (equation ^). As explained in Section 2.1, although the 
shape of the real space power spectrum in this analysis is 
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Figure 9. Likelihood contours for the two subsamples of PSCz0.6 
(galaxies with 60fi flux about 0.6Jy) split by luminosity. Here we 
have used a restricted set of modes, k max = 0.08, and a boundary 
at 200 h~ 1 Mpc. The high and low luminosity subsamples are 
created as described in the text. Contours for the low luminosity 
sub-sample are to the bottom right of the plot. The amplitude, 
Arj.i, is seen to be significantly higher for the low luminosity sub- 
sample. 
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Figure 10. Likelihood contours for the two subsamples of 
PSCz0.75 (galaxies with 60/x flux about 0.75Jy) split by luminos- 
ity. Again, we have used a restricted set of modes, k max = 0.08, 
and a boundary at 200 h^ 1 Mpc. It can be seen that by excluding 
the galaxies in PSCz with 60fim flux between 0.6Jy and 0.75Jy 
before splitting the sample by luminosity we obtain parameter 
estimates for the high/low luminosity subsamples which are con- 
sistent. 

fixed, the amplitude at a wavenumber k = 0.1 /iMpc -1 is the 
second free parameter in the likelihood fit. A wavenumber- 
dependent weighting function (equation |l3| ) has been em- 
ployed. Again, we have used a CDM-like power spectrum 
with r = 0.2 for P w (k). The normalization of P w (k) is such 
that the amplitude Aq.i is 0.18. In this minimum variance 



1 .5 - 




0.0 0.2 0.4 0.6 0.8 1.0 

Figure 11. Likelihood contours for the two parameters ft = 
f2[j- 6 /6 (where b is the IRAS bias parameter) and the amplitude 
of the power spectrum (parametrised by A 2 (fc) = k s P(k)/ (2tt 2 ) 
at a wavenumber k = 0.1 h Mpc -1 ) for the PSCz survey flux- 
limited to 0.75Jy. The outer boundary is 200 h -1 Mpc, and all 
wavenumbers up to k = 0.l3h Mpc -1 have been analysed. A 
conservative mask with extinction Ag < 0.75 has been applied, 
leaving a total of 7042 galaxies in the sample. 

weighting scheme, the exact form and amplitude of the func- 
tion P w (k) will not be critical. 

We have employed a value of a v = 420 kms -1 for the 
three-dimensional velocity dispersion of galaxies in order 
to calculate the scattering matrices and correct for non- 
linearities (equation |52| ). This is the value of a v obtained 
by adopting a pairwise velocity dispersion of 340 kms -1 
(Huchra et al. 1983, Fisher et al. 1994) and assuming that 
the small scale peculiar velocity field is incoherent. 

Figure [n] shows the contours of likelihood in the (3 — 
Ao.i plane, plotted at intervals S In/! = 0.5 for the PSCz 
survey. The recovered values of /3 and Ao.i arc: 

= 0.58 ±0.26 
Ao.i = 0.42 ±0.03. 

We have analysed PSCz0.6 to a wavenumber limit of 0.12 h 
Mpc" 1 ; this gives f3 = 0.57 ± 0.25 but A .i = 0.48 ± 0.04. 
However, as mentioned earlier, there are inconsistencies in 
detail which lead us not to trust the results from the full 0.6 
Jy sample. 

5.1 The shape and amplitude of the real— space 
power spectrum of IRAS galaxies 

The stepwise power spectrum was fitted as described in Sec- 
tion 2.2 and has been applied to the 1.2Jy survey (BHT). 
Unlike the BHT analysis, we have chosen that the right hand 
bin should begin at fc max , to show how information beyond 
the observed modes can be analysed. There is a constraint 
on this bin due to the effect of mixing, but this is quite weak 
as can be seen by the large error bar on the final point in 
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k/h Mpc -1 



Figure 12. The real-space power spectrum of the PSCz redshift 
survey in dimensionless units. The curve is a CDM model with 
T = 0.2. The estimated redshift distortion is /3 = 0.47 ± 0.16 
(conditional error). 



Figure Again, reassuringly, the redshift distortion result 
is similar to that of the two-parameter fit with 



/3 = 0.47 ±0.16, 



(24) 



although note that this is a conditional error, as are the er- 
ror bars on the power spectrum. The analysis has r max = 
200 /j _1 Mpc and fc max = 0.13/i Mpc -1 . Reducing the 
wavenumber limit to 0.12 gives a consistent result of (5 = 
0.60 ±0.18. 

The power spectrum again fits low V CDM models 
very well, although there is some evidence that the break 
is slightly sharper as suggested by Peacock (1997). It is per- 
haps possible that the final, high-A: point is suppressed by 
nonlinearities, although the random smoothing model was 
used. A particular point to note is the very tight lower limit 
on the first two bins, close to the value expected in cur- 
rent models. This demonstrates that the sample we finally 
use after the cuts discussed in Section 4.3 is uniform (in 
the sense of rms variations) to better than 3% on scales of 
500 h^Mpc. 



6 DISCUSSION AND CONCLUSIONS 

We have presented the results of a likelihood analysis of a 
spherical harmonic expansion of the PSCz galaxy survey. 
The aim has been to establish both the real-space power 
spectrum and the degree of redshift distortion from the dis- 
tribution of galaxies in redshift space. The main assumptions 
in the analysis are that structure grows via gravitational in- 
stability, and that the galaxies are linearly biased tracers 
of the underlying mass distribution, at least on the scales 
which we analyse (wavelengths > 4&h~ x Mpc). There are 
many attractive features about this analysis, and it is worth 



reiterating them here. Firstly, the expansion in spherical po- 
lar coordinates is ideal for many redshift surveys, but partic- 
ularly for the IRAS galaxy surveys, which cover a large frac- 
tion of the sky. The effects of redshift distortion are radial, 
so are naturally accounted for in this system; in addition the 
masked regions of sky, due to the galactic plane, unobserved 
regions of sky, and patches badly affected by cirrus, are spec- 
ified in terms of the angular coordinates alone. These facts 
allow a relatively simple translation from real- to redshift- 
space expansion coefficients in linear theory. This is a ma- 
jor advantage for relatively shallow IRAS surveys, since the 
linear modes typically subtend large angular scales. A gen- 
eral feature of redshift distortion studies to obtain f3 is that 
all the information is obtained from the galaxy positions, 
rather than relying on peculiar velocity measurements. The 
particular expansion we have used here is effectively a power 
spectrum analysis in spherical coordinates. Because of the 
survey geometry and redshift distortions we do not mea- 
sure the power spectrum directly, but the window functions 
converting real-space power to redshift-space expansion co- 
efficients are narrow, and we are able to measure with some 
accuracy the power spectrum on large scales between 40 and 
270 Mpc. In addition constraints on the power can be 
put on scales up to ~ 500 /i~ 1 Mpc. 

The analysis we have done follows closely the work of 
Heavens & Taylor (1995) and Ballinger, Heavens & Tay- 
lor (1995), although we have made some technical improve- 
ments in the treatment of the mask and the local group mo- 
tion, but the effects of the approximations in the previous 
papers are in fact negligibly small. 

The results we present are conservative in a num- 
ber of ways. We have firstly applied a rather conservative 
mask, and consider only regions of sky where the extinction 
Ab < 0.75. Secondly, we found that we were able to obtain 
robust results only for a subset of the catalogue cut at a 
higher flux than the formal 0.6 Jy limit. The most striking 
example of the inconsistency near the flux limit is illustrated 
in Figure [], where we have taken care to analyse two sub- 
samples with the same radial and angular distribution to 
minimise other sources of error. This graph shows why we 
do not trust the catalogue at the lowest flux limits, although 
we have not been able to pin down the source of the discrep- 
ancies. What we do find is that the results are quite stable if 
we include galaxies down to a flux limit of 0.75 Jy. The price 
we pay for these conservative choices is that the galaxy sam- 
ple is cut back to about 7000 galaxies, and the error bars 
are larger than we would have liked. However, we believe 
that the resulting estimates and errors can be trusted. In 
the first analysis, we assume that the power spectrum has 
the shape of a T = 0.2 CDM spectrum, at least on scales 
larger than 40 h~ Mpc, and estimate simultaneously (3 and 
the amplitude (characterised by the dimensionless power at 
0.1 /i Mpc). This analysis yields 



(3 = 0.58 ±0.26 
Ao.i = 0.42 ±0.03. 



(25) 



If we allow the power spectrum to be free (or at least charac- 
terized by a large number of parameters), the results change 
only slightly: 



0.47 ±0.16 
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Ao.i = 0.47 ±0.03 

(26) 

but note that in this case the quoted errors are conditional, 
and are not reflective of the true marginal errors. 

The value of which we derive is lower than the values 
of 1.1 and 1.0 obtained by similar analysis of the IRAS 1.2 
Jy catalogue (Heavens & Taylor 1995, Ballinger, Heavens 
& Taylor 1995), and by the analysis, also using spherical 
harmonics, of Fisher, Scharf & Lahav (1994). However, the 
errors in our previous analyses are large enough to allow f3 to 
be as low as 0.5. Our result is consistent with determinations 
based on the plane-parallel approximation (e.g. Cole, Fisher 
& Weinberg 1994, Cole, Fisher & Weinberg 1995) and corre- 
lation function studies (Hamilton 1992). A full discussion of 
results in this area can be found in Hamilton (1997). Meth- 
ods based on comparison between velocity fields and IRAS 
density fields also show a (formally inconsistent) spread of 
results, from 0.49 ± 0.07 (Willick et al. 1997) to 0.89 ± 0.12 
(Sigad et al. 1998). One possible way to reconcile discrepant 
results in general is if bias is scale-dependent (e.g. Coles 
1993, Mann, Peacock & Heavens 1998), and the measure- 
ments are effectively probing different scales, but it is not 
clear whether this is appropriate here. 

Since we compute the power spectrum and the red- 
shift distortion simultaneously, and assuming the linear bias 
model, we are able to express the answer in terms of the 
strength of clustering of the mass, as has been done from 
the abundance of clusters (White, Efstathiou & Frenk 1993). 
We express the amplitude of matter clustering in terms of 
the dimensionless matter power spectrum at a wavenumber 
of 0.1 h Mpc" 1 : 



Ao.i, p = 0.24fio 6 (27) 



with an error of about 50%. If we convert this to the r.m.s. 
fluctuations in an 8 /i~ 1 Mpc sphere, assuming a Y — 0.2 
shape parameter (using the approximation equation jl6]) , we 
get o"8 = 0.38f2,7°' 6 , again with an error of about 50%. This 
rises to 0.57 if V = 0.5. 

What should we conclude from this? The simplest, al- 
though not the only, interpretation is that IRAS galaxies 
are unbiased, and fto — 0.4. The strongest indications that 
the Universe might be Einstein-de Sitter have in the past 
come from the high values of /3 from peculiar velocity and 
redshift distortion studies. The picture is still not clear, but 
the evidence is shifting away from a high f3 towards low. In 
any case, it is impossible to make firm statements about Qq 
based on these studies because of our ignorance of the bias 
parameter, and we shall almost certainly have to wait until 
large redshift surveys of 10 5 - 10 6 galaxies are completed, 
when higher-order statistics will lift the degeneracy between 
Qo and b (e.g. Matarrese, Verde & Heavens 1997, Verde et al. 
1998). Indications of a low value of Slo come from the super- 
nova Hubble diagram (Perlmutter et al. 1998, Riess et al. 
1998) combined with rather uncertain flatness constraints 
from the microwave background, although the prospects are 
good (Tegmark, Eisenstein & Hu 1998). 
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APPENDIX A: OUTLINE OF THE SPHERICAL HARMONIC FORMALISM 

In this Appendix we outline the spherical harmonic formalism used to estimate the power spectrum and redshift space 
distortion in the PSCz. The general derivation is similar to that given in HT, but is somewhat more comprehensive. 

A.l Definitions 

Following HT we expand the density field of the galaxy distribution in spherical harmonics, Y( m , and spherical Bessel functions, 

je, 

Pimn = c ln J d 3 s p(s)w(s)j e {k ln s) Y e * m (9, 4>) , (28) 

where w(s) is an adjustable weighting function and s is the redshift-space position variable. ctn are normalisation constants, 
and k£ n are discrete wavenumbers (see A2). The inverse transformation is 

p(s) = ^ c£ n pi mn je(ke„s)Y em (0, 0). (29) 

imn 

We use the Binney and Tremaine (1987) definition of spherical harmonics which uses the Condon-Short ley phase condition; 



/ 21+1 (I - \rn\y. D , m |, m4l J (-l) m m>0 

where Pj m '(x) is the Legendre polynomial defined (for positive m only) by 

im 

pr( I ) = (i- I r /2 ^ftw. (31) 

Note that some definitions add an extra ( — l) m term (e.g. Press et al. 1980). With the Binney and Tremaine definition the 
Hermitian relation, 

Yf m = (-l) m y/ ra *, (32) 
is automatically satisfied. 

A. 2 Normalisation and boundary conditions 

The ct n are normalization constants introduced to preserve the normalisation of the spherical Bessel function in a finite 
sphere. They are derived via the application of the orthogonality relations for spherical Bessel functions, for given boundary 
conditions (see HT). The boundary conditions used for this analysis are Neumann boundary conditions (i.e. that the first 
derivative of the potential field normal to the boundary is zero). This implies that the velocity field is zero on the boundary i.e. 
the boundary is undistorted. We therefore need not include surface effects. The boundary condition gives us the wavenumbers 

^-ji (k en r) | r =r max = 0, (33) 
ar 

and the normalization constants 

kin /oa\ 

Ctn = , (34) 



\je(ktnr m &x)\ 



4 ' In' max 



A. 3 Redshift Distortions 

As described in the main text, the transformation to redshift space can be carried out by the identity d 3 s p(s) — d?r p(r), 
leaving redshifted coordinates in the weighting and spherical Bessel functions. HT then expand w(s)je(ki„s) in equation ( |28[ ) 
to first order in s — r, 

w{s)ji{k ln s) ~ w(r)j e (k tn r) + u(r)-^- [w(r)j e (k en r)] , (35) 

ar 

where u(r) is the radial component of the peculiar velocity field, v(r) — ulg- Linear theory provides a relation between the 
peculiar velocity field v(r) and the density field via the continuity equation 

V.v(r) = ~H Q ^ s 8 p (r) (36) 
to first order. In Appendix A of HT it is shown that the radial component of this field is 
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u(r) = Qq 6 \_. cink e ^ 5tmnj't{ktnr)Yt m (6, <f>) - dlg ■ r, (37) 



tran 

where j'(z) = dje(z)/dz and 8e mn is the transform of the galaxy number over-density field 8(r). Note that in this expression we 
have added the contribution from the Local Group dipole. This was neglected in the analysis of HT since it only contributed 
significantly to the 1 = 1 mode, which was not included in the analysis. We include it here for completeness. 
Substituting equation ( pTj ) into equation ( |35| ) gives 

w(s)j e (k en s) ~ w(r)j e (k en r) + ^ Qwfciw&Wn'ji'Cfovr)^ [ w ( r )je(ke n r)]Y e > m ' - (i>lg ■ r)^ [w(r)j>(fc fn r)] . (38) 

A. 4 Transforming real space harmonic modes to redshift space 

From the definition of perturbations in the galaxy density field we have 

p(r) = po(r) [1 + S(r)] = p (r) I 1 + c ln 8i> mrl j t (k e „r) Y lm (9, 0) , (39) 



where po{r) is the observed mean density of galaxies in the survey. 
Substituting equations (B8I) and (133) into equation (E8h we find 



Ptmn = {po) tmn + ^ ( $ «W + W™') 8$L'- (40) 



Vm'n 1 

The mean-field harmonics are defined as 

(Po)imn = y d\ p {r)w(r)ji {k (n r) Y e * m (9, <f>) - c;„ / d 3 r («lg ■ r)po(r)j^; [iu(r)j* {k tn r)\ Y/ m (9, 4>) , (41) 

where we have incorporated the the contribution of the dipole term, since this is independent of 5(r). In the rest of the 
analysis we assume a Local Group velocityof 622 km s _1 towards t — 277°, b = 30°. 

The remaining matrices in equation (HOI) account for the mixing of modes due to the angular and radial geometry of the 
survey: 

*«W = ctnCe' n > Jd 3 r po(r)w(r)j e (k ln r) Y e * m (9, </>) j e , (fo/ n /r) Y t i m , {9, </>) (42) 
and the angular geometry and radial distortion: 

V«W = Ce "^ e ' n ' / d 3 rp (r)-^- [w(r)je (k in r)] Y e * m (9, (j>) j' e , (k e , n ,r) Y e , m , (9, </>) . (43) 



kj, n , J dr 

Defining the observable data vector De mn — (pemn — (po)tmn)/p, and <5, the vector of underlying density field values, 
equation mm can be written in a compact form as 

£» = (* + pV)8. (44) 
where we have adopted a matrix notation. 

A. 5 Separation into angular and radial terms 

The mean observed density field, po(i*), may be separated into angular and radial components 

po(T)=M(6,4>)po(r), (45) 

with M(6, (j>) being the sky mask of the survey and po(r) containing the radial selection function. The 'transition matrices' <& 
and V now separate to become 

®ee'nn> = W U i (46) 

U'nn' — W W V W I 4 'J 

where 

=ct n ct>n> / r 2 drp (r)w(r)j e ,(k e ,„,r)j e (k en r), (48) 
Jo 

V" n > = / r 2 drp (r)^- [w{r)j l (k ln r)] ^-j t , (k e , n ,r) (49) 

L'fi' JO ' 

and 
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w%r = / y £ t (n)M{n)Y; m {n)dQ. (50) 



A. 6 "Fingers of god" and nonlinearities 

Throughout the paper we have also applied the correction (described by HT) for the effect of the small scale non-linear 
peculiar velocity field (see their Appendix B). In this procedure, it is assumed that the real-space position of a galaxy is 
perturbed slightly by a peculiar velocity such that the radial redshift-space position becomes 

s(r) — » s'(r) = s(r) + e(r). (51) 

The quantity e(r) is a random variable drawn from a Maxwellian distribution with mean zero and (e(r)e(r')) = 
(a^/3HQ)6D(r — r') and where a v is the 3-d velocity dispersion. Errors in the measurement of redshift, cr 2 . can be incor- 
porated into the formalism by adding the variances in quadrature; crf o4o ; = (ff^/3) + o\. Since this distortion is incoherent, 
the effect is to convolve each mode with a radial Gaussian distribution. From this HT derived a scattering matrix S defined 
by 

(r-v) 2 
V5" 



' ' -je(k en r)ji> (kt> n >y) rdrydy (52) 



Vn J J V27ra 

which is to be convolved with the density field. S K is the Kronecker delta function. The amplitudes, Di mn , corrected for a 
small-scale incoherent velocity field are now given by 

D lmn = S«™„ n //WP™ +pvjtn? ) iftnw, (53) 

or in matrix notation, 

D = SW(® + (3V)5. (54) 
A. 7 Correlations 

HT assumed that the sky mask was azimuthally symmetric, making the mask mixing matrices W real. In the current analysis, 
we drop this assumption as the sky mask for the PSCz survey is not symmetric (see Figure |?j). This complicates the analysis, 
since the noise matrix is no longer real. We can avoid this problem by constructing the covariance matrix for the real and 
imaginary parts of D. Since these are also Gaussian random, there are now four terms to the covariance matrix C , 

CW = {OkDOte-D*), (55) 

C»9m = (ZhDQmD*), (56) 

CWfc = {SmDVhD*}, (57) 

= (QmDQmZ?'), (58) 

where 5Rfe signifies the real part and S»n the imaginary part. These terms are given by the following summations: 

C'QmGm = Csfegfc, (60) 

C*s*n = i J] [5fe(<£>£ + /3V M a )an(C + PV?) - 3m(<^ + + /?K*)] P(*„) + N^, (61) 

a 

C'omSKte = — Csteam, (62) 
where /j = (£mn), v — (i'm'n') and a = (£"m"n"). The noise matrix, iV, is given by the equation 

N„ v =c en c e , n , / p (r)w 2 (r)j f {k tn r)j t > (k e , n ,r)r 2 dr / P M [^(ft)] M{Q)V V [Y? m , (fi)] dfi (63) 



as in HT, where the Vy,, v represent the real or imaginary parts of the spherical harmonics according to whether we are 
using the real or imaginary parts of D^. v . In Appendix B we discuss in more detail the construction of the mask and its 
decomposition into correlations between real and imaginary parts of the field D. 

Now we have the covariance matrix, a likelihood functional is constructed and we use an assumed real-space power 
spectrum P(k) in equations (jH^) to (^). We maximize the likelihood for the value of fi and the amplitude of the power 
spectrum, as described in the main text. The amplitude of the power spectrum is defined in terms of Ao.i where A 2 (fc) = 
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fc 3 P(k)/2ir 2 . The likelihood function (unnormalised) we use is 
£[D\f3,P{k ln )] = [det(C ste!R! )det(C S K tam )det(Ccj m!R! )det(C 3in c Jm )]-i 
x exp (-hfcD'C^-SkD - ^D*C£riJkD) 

x exp^i^C^QmD-iomD'C^QmD) . (64) 



APPENDIX B: CONSTRUCTION OF THE ANGULAR WINDOW MATRIX 

In this appendix we describe in detail the formulation and calculation of the angular window matrix. Since we drop the 
assumption that the angular mask is axisymmetric, calculation of the data covariance matrices is somewhat more complicated 
than that found in HT (see Appendix A. 7). In this Appendix we define the window transition matrix, W , and describe how 
we can calculate it efficiently with the use of Clebsch-Gordan, or Wigner 3j coefficients. We then show that formally the data 
covariance matrix is imaginary, due to the non-axisymmetric mask. We expand the window matrix and show how the data 
covariance matrix can be made positive definite by taking the covariances of the real and imaginary parts of the data vector, 
D. Finally, we present the explicit expressions used to calculate the full covariance matrix, including the noise terms. 



B.l Definitions and conventions 

An angular mask function is defined as 

,._J I if observable area of sky, , . 

I otherwise 

From this we can define the transition matrix 

w%r' = f Y e r'(n)M(n)Y e m *{n)dn, (66) 

J Air 

In this Section we shall abbreviate our notation by dropping I, m subscripts and J dQ except where needed for clarity. Note 
our convention for dashes and complex conjugates in the definition of W . The W matrix operates on the density field by 
mixing modes; 

D'tmn = ^ Ww™ D t i m , n . (67) 
I'm' 

Rather than directly performing the integration in equation (|6(]), which is a prohibitively slow process for large £, it is faster to 
decompose the calculation into a two-step process, as outlined by Scharf et al. (1992). First one calculates the single transform 
of the mask function 



Mi m = / M{Q)Y e m *{Q)dQ, (68) 
which can then be use to calculate the full transition matrix by the transformation (Y' = Y™ etc.) 



W%r = 2^ MP, I Y'Y"Y* . (69) 

t"m" 



B.2 Clebsch— Gordan & Wigner 3j coefficients 

The angular integral, J Y'Y"Y* can be calculated analytically by the Clebsch-Gordan or Wigner 3j matrices (Edmonds 
1957). The Wigner 3j symbols are related to the Clebsch-Gordan coefficients by 

( £ £ ', l "„ \ = (-l) e - e '- m "(2£" + l)- 1/2 (£m£'m'\££'£"-m") (70) 
\ m m m I 

Here we choose to express relations in terms of the 3j matrices. Hence 

fw'v" _ / (2l + l)(2l' + l)(2l" + l) (££'£" \ ( £ £' £" \ 
J YY Y "V 4^ ^ J { m m' m" J 

with the triangle conditions m + m' + m" — and \£ — £'\ < I" < (£ + £'), plus permutations. Integrals containing complex 
conjugates of Y's can be derived using the hermicity of the spherical harmonics (equation (^)). 
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The Wigner 3j symbols are defined by Landau & Lifshitz (1977) 



771 m' m" J °mi+m 2 ,m { J- J ^ + £> + £» + 

[(£ + m)\{£ - m)\{£' + m')\{£' - m')\(l" + m")\{l" + m")}] 1/2 
^2(-l) z [z](e +£'- I" - z)l{£ - to - z)\{£' + to' - z)\ 

Z 

{£" -£' + m + z)\{£" -£' - m + z)!]' 1 (71) 
where the summation over z is zero when any factorial terms in the summation are negative. When m = m' — m" — then 
I £' £" \ L/3 I (L - 2£)\(L - 2£')\(L - 2£")\ (L/2)\ 



t" \ , .si/a l (L-2£)\{L-2£')\{L-2£")\ 
) V 



J V (L + lV- {L/2-£)\{L/2-£')\{L/2-£")V 

I odd. (72) 

£ £ \ 

where L = £ + £' + £" . Practically it is easier to calculate these terms from [ I using the recursive relations 



(73) 



£ f f \ _ (L-2£' -l)(L-2£" + 2) ( £ £' + 1 £" - 1 
J^Y (L-2£')(i-2^" + l) V 

With these relations, the computation time for the mask mixing matrix is reduced from days to a few minutes for £ up to 100. 
B.3 Covariance Matrix Terms 

Armed with the basic definitions of the window function transformation, and methods to generate it efficiently, we now begin 
to enumerate the various terms that will appear in the data covariance matrix. As stated above the problem of constructing 
a data covariance matrix is complicated by the combined presence of noise and a non-axisymmetric mask, leading to an 
imaginary covariance matrix. This can be dealt with in linear theory by decomposing the harmonic modes into real and 
imaginary parts and calculating their covariances. We begin with a formal derivation of the data covariance matrix. 
The covariance matrix of the observed modes can be written in matrix notation as (c.f. equation (@)) 

C = (DD 1 ) = SW{® + j3V)P{® + PVY W^S* + RW, (74) 

where P = dia,g\P(kz„)] and the shot noise term is given by N = RW , which we have decomposed into radial and angular 
parts (c.f. equation (p3|)), 



R = ce n Ci' n ' I po(r)w {r)ji (ki„r) j e , (k e , n ,r) r dr (75) 
and 

W = I Y'MY*. (76) 



Note that here we use the full window transition matrix. 

Expanding W into real and imaginary parts, W = $hW + i^smW , yields 

C = S?fhW($ + f3V)P{® + pV) t $kW t S t + S3mW(& + ,0V)P(* + /3V) t '^mW t S t + R$hW 

+ iS^W(® + /3V)P($ + l3V) t m:W t S t -iS^W(® + l3V)P(® + l3V) t QmW t S t + iR^W (77) 

where the real and imaginary parts of W can also be written 

VhW = SymW, (78) 
and 

3mW = -{-AsymW, (79) 
i 

where the operators Sym and Asym symmetrize and anti-symmetrize the matrices. 



B.4 Expansion of the window matrix 

Although imaginary terms exist in the data covariance matrix, we can use the property of a Gaussian field that the real 
and imaginary parts will be independent and each is independently Gaussian distributed to construct total the probability 
distribution function. We later show that the covariance matrix can be further expanded into independent terms that will 
allow us to define a total covariance matrix which is both real and symmetric. 
Expanding the Y's into real and imaginary parts we find that 
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(80) 



where we introduce the notation WShW = J ShYMShY', SmSimW = f <3mYM%mY', JfcQmW = / 3hYM%mY' and 
•vmSkW = / QmYMVkY'. The real and imaginary parts of the spherical harmonics can be found from 



and 



9my = — (Y 



Y* 



(81) 



(82) 



The terms SffeSfSsW can then be calculated from the full window matrix by 



r 
i 

47 



YMY' + I Y* MY' + / YMY" + / Y* MY'* 



Y* MY' 
YMY' - 



YMY'* 



YMY' 



YMY'* + / Y* MY' 



Y* MY 1 * 
Y* MY'* 



(83) 



where the various permutations can be calculated from the window function matrix using the Hermicity relations. These are 
the terms that will appear in the noise matrix, equation (^). Note the ordering of dashed terms as these terms are not in 
general symmetric in their indices. Clearly then 3hW = 3ffe5SsW^ + ^Sni^smW and Ssm W = 3hQmW — ^sn^ShW . 



B.3 Transpose relations 

Here we list without proof the transpose relations for the three window matrices. 

Sfe3feW = {?MkW) T , 
Qm^smW = (Qm?smW) T , 

Qm'SkW = (3fe3mW) T . (84) 

Given the symmetry relation Qm?ShW = (3hQmW) T , we only need explicitly to calculate the three window matrices 3h?lhW , 
Qm^smW and Qm^hW . 



B.4 Parity Relations 

Here we list without proof the parity relations (m 



m. — m 



= (-l) m $h%eW„ 



m^ — m 



m. — m 



= (-ir 



= (-l) m ' +1 9mQmW„ 
= (-l) m+1 Qm3mW m 



(-1) 



m + m CVnCi 



QmQmWrr 



-m) for the three basic window matrices. 



(85) 



(86) 



m. — m' 



(-l) m 9hi3feW mTr . 
(-l) m+1 9m3feM/„ 



(87) 



B.5 Expansion of covariance matrix 

We can re-arrange the covariance function by reducing the data vector into real and imaginary parts; D — 5fife.D + i^mD. 
This is convenient because our hypothesis is that D is a Gaussian random variable and so the terms IShD and ^sraD will also 
be Gaussian random variables. In this case we need terms like {$hD $hD'). This leads us to 

C = (QteD $hD n ) + (SmD SmD n }) + i((3mD $hD tt ) - (3kD <3mD n )). (88) 
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We can further expand terms by expanding the real and imaginary terms of the W matrix and the underlying field, 5. Hence 
$hD — $hW $h5 — QmW ^smS and QmZ) = QmW $h8 + $hW QmS. Under the assumption of statistical independence between 
the real and imaginary <5's we find 

(MsZ? tfteD 1 ) = ^S^W($ + pV)P($ + (3V) t ^W t S t + S^W($ + /3V)P($ + i3Vy^W t S t + R^$hW, 
(SmDZmD*) = $hW($ + /3V)P{<f> + (3V) f WtW* S l + S <3mW($ + f3V)P{® + j3Vf ^nW'S 1 + R QmQmW , 
(XkDSmD 1 ) = ^S$kW(& + /3V)P{& + fiV)* &nW t S t - S3mW{& + /3V)P(# + /3V)* SteW'S* - RShQmW , 
{SmDZteD*) = ^S3mW(& + /3V)P(& + /3V)* OfeW'S* - S$kW{& + /3V)P(* + /3V)' OmW'S* - RSmVhW. 

(89) 

Substituting these expressions into equation Q88p we reproduce equation (Q). In addition we have shown that each of these 
terms are independent and Gaussian distributed. Therefore we can construct a real covariance matrix; 

C=( Z**' ) (90) 

From the symmetry relations given earlier we see that this is a real, symmetric matrix. 
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